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We propose a new algorithm of the finite lattice method to generate the high-temperature series 
for the Ising model in three dimensions. It enables us to extend the series for the free energy of the 
simple cubic lattice from the previous series of 26th order to 46th order in the inverse temperature. 
The obtained series give the estimate of the critical exponent for the specific heat in high precision. 
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The finite lattice method [y, |^ is a powerful tool to 
generate the high- and low- temperature series for the 
spin models in the infinite volume limit. It avoids the 
tedious work of counting all the diagrams in the graph- 
ical method and reduce the problem to the calculation 
of the partition function. In two dimensions the total 
amount of the calculations for the finite lattice method 
increases exponentially with the maximum order N of 
the series. On the other hand in three dimensions the to- 
tal amount of the calculations grows exponentially with 
N"^ and except for some cases § |, |, |, | ||, |l^, |ll|, |l| 
many of the expansion series have been calculated by 
the graphical method. Here we present a new algorithm 
of the finite lattice method for the high-temperature ex- 
pansion in three dimensions in which the total amount 
of the calculation increase approximately exponentially 
with TV log N and that enables us to generate the series 
to much higher orders than not only the standard algo- 
rithm of the finite lattice method but also the graphical 
method. 

In the finite lattice method to generate the high- 
temperature series for the free energy in three dimensions 
we calculate the partition function Z{lx x ly x Iz) for the 
finite size lattices with 2{lx -\- ly Iz) < N and define 
recursively Q 



(l){lx Xly X Iz) = log [Z{lx XlyX Iz)] 
C</x,i;</y,C<'.: 



(1) 



Here we use the notation for the lattice size such that the 
1x1x1 lattice means the unit cube composed of 2 x 2 x 2 
sites. The Boltzmann factor for each bond connecting the 
nearest neighbor sites k and k' is expressed as 



exp {f3skSk') = cosh(/3) (1 -f tskSk') , 



(2) 



with t — tanh (/?) . We define the bond configuration as 
the set of bonds to which the factor tskSk' in (||) is as- 
signed while the factor 1 is assigned to the other bonds 



of the finite size lattice. Non-vanishing contribution to 
the partition function comes only from the bond con- 
figuration in which the bonds form one or more closed 
loops. Each of the closed loops is a polymer in the stan- 
dard cluster expansion Ip^ . Then the Taylor expansion 
of X ly X Iz) with respect to t includes the contri- 
bution from all the clusters of polymers in the standard 
cluster expansion that can be embedded into the lattice 
of Ix X ly X Iz but cannot be embedded into any of its 
rectangular sub-lattices l'^ xl'yX l'^ ||] . The expansion se- 
ries of the free energy density in the infinite volume limit 
is given by 



/ = E "^('^ xly X Iz) 

2{l^+ly+h)<N 



(3) 



The expansion series of (f>{lx x ly x Iz) starts from the 
term of with n — 2(1^ + ly + Iz), which comes from 
the cluster of a single polymer (one closed loop of bonds) 
that have two intersections with any plane perpendicular 
to the lattice bonds. Thus it is enough to restrict the 
lattice sizes for the summation in (^) to those that satisfy 
2{lx + ly + Iz) < N to obtain the series for / to order . 

In the standard algorithm of the finite lattice method 
the full partition function for the finite size lattice is cal- 
culated with all the bond configurations taken into ac- 
count. The point of the new algorithm is that, in order 
to obtain the series to a given order, however, it is enough 
to consider only a restricted number of bond configura- 
tions. Let us consider the anisotropic model of the simple 
cubic Ising model with (3i = Ji/ksT and ti = tanh(/3i) 
{i = x,y,z). To obtain the series for (j)(lx x ly x Iz) to 
order Nz — 2lz + ANz in tz we introduce in the new 
algorithm ip{lx x ly x lz,ANz) defined recursively by 

(t>{lx XlyX lz,ANz) - \og[Z{lx XlyX lz,ANz)] 

J2 mxl'yXl',,ANz). (4) 

l'^<l^,l'y<ly:l'„<lz, 

I'^+l'y+l'zT'^l^+ly+l^ 



Here the partition function Z{lx x ly x lz,ANz) is cal- 
culated only with the bond configurations taken into ac- 
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(a) (b) 




FIG. 1: Examples of the bond configuration for Z(lx x ly x 
l^,AN^) witli Ix = 3,ly ^ 3,lz = 4, AiVz = 2. Example (a) 
has Wzi = 2, nz4 — 2, n^a — 4, and = 2 and it should be 
taken into account, while (b) has Uzi — 0, nz4 — 4, n^a = 4, 
and = 2 and it should be neglected. 

count that have orders Uzi in for the i-th layer per- 
pendicular to the z-axis (i — 1,2, ■ ■ ■ ,1^) satisfying 

2max(n^„2) < 2/^ + AiV^ . (5) 

Examples of the bond configuration for Z{lx x ly x 
Iz, ^Nz) that should be taken into account or should be 
neglected can be seen in Fig. |l|. 

It is easy to prove that any bond configuration for the 
partition function that has X^iLi "^zi £ 2/^ + ATV^ but 
Yl\'Li max(nzi, 2) > 2/z + AA^z does not contribute to the 
left hand side of in the order lower than or equal to 
Nz = 21 z + AA^z. For such a bond configuration at least 
one of the nziS should be zero, so they are disconnected 
configuration(composed of more than one polymer) or 
they can be embedded into a rectangular sub-lattice of 
the Ix X ly X Iz lattice and in either case they do not 
contribute to (l){lx xlyXlz, ^Nz) in the order lower than 
or equal to Nz = 21 z -\- AA^^. They can contribute to 
4>{lx X ly X lz,ANz) in higher order than Nz = 21 z + 
ANz by constituting the connected cluster of polymers 
together with other configurations with n' zi > 2 for the 
layer i with Uzi = . Then total order of the cluster of 
polymers has the order higher than Nz = 21 z + ANz ■ 

The contribution of the bond configuration with {uzi} 
to the partition function of the finite size lattice can be 
calculated by the transfer matrix formalism as 

Z{{nz^}) = V^.nil^'Vn.n^l'^ ' ' ■^""^^.,.,o • (6) 

Here Vj.j.^-^ is the transfer matrix element with incom- 
ing Uzi spins and outgoing Uz i+i spins and the summa- 
tions over the spin locations ji, j2, • • • of the nzi,nz2, ■ ' ' 
spins, respectively, are assumed in the right hand side 
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FIG. 2: Total amount of the calculation to generate the 
series to order A'' by the old algorithm and the new algorithm 
of the finite lattice method. 



of (|6|). This transfer matrix element itself is the parti- 
tion function in two dimensions with Uzi + Uz i+i spins 
attached, which can be calculated to any order in 
and ty efhciently by the site-by-site construction!^^ p^ . 
The amount of the calculation for each transfer ma- 
trix element is proportional to the combinatorial factor 
C( (Ix + l){ly + 1) , rizi + Uz i+i), 2'- and IJy, which are 
the number of the cases for attaching the Uzi -f Uz i+i 
spins to the [Ix + l){ly + 1) sites, the number of states 
in site-by-site construction for the partition function of 
the Ising model in two dimensions and the number of the 
relevant bonds, respectively. 

To obtain the series to order A^ in the isotropic model, 
we calculate the expansion series for each of the (f){lx x 
ly X IzYs defined by (|l]) in the anisotropic model to or- 
der t^-tyH^- with Nx + Ny + Nz = N using the new 
algorithm described above and set t^, = iy = = i fi- 
nally. When we use the new algorithm for each lattice 
size Ix X ly X Iz and each of the orders Nx, Ny, Nz, we 
can make the simultaneous exchange of the lattice axes 
and corresponding orders as {Ix x ly xlz', Nx, Ny, Nz) 
(Ix xlzX ly;Nx,Nz,Ny) or {Ix x ly x h; Nx, Ny, Nz) 
{ly X Iz X Ix; Ny, Nz, Nx) that have the chance to reduce 
the amount of the calculation for the transfer matrix el- 
ements. 

We estimate the total amount of the calculation time 
Tcalc (A^) to generate the free energy series to order A^ 
by listing up all the transfer matrices needed and sum- 
ming up the estimated time to calculate each of the ma- 
trices, which is plotted in Fig. || together with the es- 
timated total calculation time for the old algorithm of 
the finite lattice method. The numerical estimation for 
the new algorithm agrees well with the actually used 
calculation time for A^ < 46. We see that the calcu- 
lation time for the old algorithm grows exponentially 
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TABLE I: High-temperature expansion coefficients for the 
free energy density of the simple cubic Ising model. 
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8410722262379235048686604 / 1 1 
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14120314204713719766888210 



with iV^, while that for the new algorithm can be fitted 
by ^ exp (BiV log TV + C7V + D log TV) with B - 0.15, 
C ^ —0.24 and D ^ 3.6. We can simply understand 
this exp(i3A''log + ...) behavior as follows. The max- 
imum size of the lattice to be taken into account is 
Ix = ly ~ Iz = N/6, for which AA^^ = and we have 
only to consider the bond configurations with n^i — 2 
for all z = 1, 2, • • • , The largest amount of the calcu- 
lation is to be paid for the partition function of the lattice 
that has smaller size of Ix ^ ly ^ Iz ^ N/12 for which 
the maximum of Uzi + nz i+i is about A^/6 and the prod- 
uct of the above factors is approximately proportional to 
exp {BN logN + CN + D log N) with 5 = 1/6 for large 
N. 

Using the new algorithm of the finite lattice method we 
have calculated the high-temperature series for the free 
energy density of the simple cubic Ising model to order 

= 46. The coefficients of the obtained series 

/ = 3 cosh + ^ a„i" ; t = tanh (/?) (7) 

n 

are listed in Table |. They agree with those given by 
Bhanot et alQ to order TV = 24 and by Guttmann and 
EntingQ to order TV = 26. We have add ten new terms to 
these previous series. To obtain the series of order TV = 
26 in our new algorithm we needed only the computer 
memory of 1 MBytes and the CPU-time of 5 minutes in 
a standard PC and to order TV = 46 we used the total 
memory of 2 GBytes and the CPU time of 25 hours in 
CP-PACS at Tsukuba University. 
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FIG. 3: Critical exponent a versus the critical value Pc 
for the inhomogeneous differential approximants of the first 
order. 
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FIG. 4: Sequence a„ plotted versus 1/n. 



Following is the result of the preliminary analysis of 
the series for the specific heat C{t) — X^n'-"^"- 
plot in Fig. 1^ the critical exponent a versus the crit- 
ical value /3c for the first order inhomogeneous differ- 
ential approximants of the series of TV = 40-46. From 
the linear dependence of a on f]^ we find a = 0.1045(1) 
or a = 0.1077(2) at the value /3c = 0.22165459(10) |l| 
or Pc = 0.2216595(15) (l^ obtained in the recent Monte 
Carlo simulations. The second order inhomogeneous dif- 
ferential approximants give almost the same value. 
We also made the ratio analysis. The sequence an = 
(icC2n/c2n-2 " l)'^' + 1, which is expected to behave as 
a + b/n^ + c/n + ■ ■ ■ with the correction-to-scaling ex- 
ponent A ~ 0.5, is plotted versus 1/n in Fig. ^ for 
/3c = 0.22165459. The sequence for n = 14-23 given by 
the new 10 coefficients has a bit different slope from the 
sequence for n < 13 given by the previously obtained 
series. Three-parameter fitting of the newly obtained 
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sequence for (3^ = 0.22165459 gives a = 0.1036(10), 
b = -0.007(10) and c = 0.17(2) for A = 0.5. As for 
the case of Pc = 0.2216595 it gives a = 0.1082(14), 
b = -0.033(10) and c = 0.21(2). These estimated val- 
ues of a by the inhomogeneous differential approxima- 
tion and by the ratio method are not inconsistent with 
each other. We see, on the other hand, that the esti- 
mated values of a depends crucially on the value of /3c 
and in order to determine a precisely in these biased 
method we need more precise value of /3c. We also notice 
that the correction-to-scaling term is very small, which 
was already pointed out in the analysis of the shorter 
series @, ||]. 

From the hyperscaling relation a — 2 ~ dv the high- 
temperature series for the second moment correlation 
length gives a = 0.1088(39)@ and a = 0.1096(5)[^ 
and the e-expansion gives a = 0.1088(39)||2^. More di- 
rect estimation of a = 0.110(2) was also given by the 
high-temperature series for the magnetic susceptibility of 
the antiferromagnetic critical point ||l^ . We note that our 
estimated value of a using /3c — 0.22165459 is not consis- 
tent with these values but the value using /3c = 0.2216595 
is rather closer to these values. 

The basic idea presented here in the new algorithm 
of the finite lattice method can be applied to the high- 
temperature expansion of other quantities such as the 
magnetic susceptibility and the correlation length for the 
Ising model in three dimensions and it can also be applied 
to the models with continuous spin variables such as the 
XY model in three dimensions. Furthermore the idea can 
be used in the low-temperature expansion for the spin 
models in three dimensions. We can expect that it will 
enable us to generate the expansion series in much higher 
orders compared with the presently available series. 
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